Facilitated spin models of dissipative quantum glasses 
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We introduce a class of dissipative quantum spin models with local interactions and without 
quenched disorder that show glassy behaviour. These models are the quantum analogs of the 
classical facilitated spin models. Just like their classical counterparts, quantum facilitated models 
display complex glassy dynamics despite the fact that their stationary state is essentially trivial. In 
these systems, dynamical arrest is a consequence of kinetic constraints and not of static ordering. 
These models display a quantum version of dynamic heterogeneity: the dynamics towards relaxation 
is spatially correlated despite the absence of static correlations. Associated dynamical fluctuation 
phenomena such as decoupling of timescales is also observed. Moreover, we find that close to the 
classical limit quantum fluctuations can enhance glassiness, as recently reported for quantum liquids. 



Introduction. A central problem in condensed-matter 
science is that of the glass transition. Many body sys- 
tems with excluded volume interactions, such as molec- 
ular liquids, experience pronounced dynamical slowdown 
at high densities and/or low temperatures, to the ex- 
tent that they eventually cease to relax and form the 
amorphous solid we call glass. The glass transition as 
observed experimentally is not a phase transition but a 
very rapid kinetic arrest. At low enough temperature or 
high enough density glass formers relax too slowly to be 
observed experimentally in equilibrium and thus behave 
as (non-equilibrium) solids. This solidification occurs in 
the absence of any evident structural ordering, in con- 
trast to more conventional condensed matter systems: in 
glass formers thermodynamics changes apparently very 
little but dynamics changes dramatically. Dynamic ar- 
rest like that of glasses is a generic phenomenon in con- 
densed matter; for recent reviews see [IHI]- 

While the first hallmark of glass formers is kinetic ar- 
rest, the second is dynamical heterogeneity jj]. Glass for- 
mers appear structurally homogeneous, but their dynam- 
ics is highly heterogeneous: as they slow down spatial dy- 
namical correlations emerge and these become more pro- 
nounced the longer the relaxation times. One theoretical 
perspective on glasses in which dynamical heterogeneity 
appears naturally is that of dynamical facilitation, which 
posits that the origin of glassy slowing down is not to be 
found in thermodynamic ordering [6 but in effective con- 
straints in the dynamics (see [5] for a review) . From this 
perspective, slowdown, heterogeneity and other fluctua- 
tion features of glasses are rooted in the complex struc- 
ture of trajectory space. This theory has emerged from 
the study of a class of idealized lattice systems, so-called 
kinetically constrained models, of which the simplest rep- 
resentatives are the facilitated spin models [7J- 

Quantum glasses, just like their classical counterparts, 
are of much current interest, among other reasons due 
to their relevance to issues like supersolidity jS], quan- 
tum annealing [S], glassiness in electronic systems [ID] . 
thermalization and many-body localization |llj . and ar- 
rest in quantum fluids [12]. Central questions in quan- 
tum many-body systems which display glassy behaviour 
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FIG. 1: Quantum facilitated models. A lattice of two-level sys- 
tems with coherent and incoherent excitation/de-excitation. 
Interaction is via kinetic constraints fi which make the rates 
on site i dependent on the state of its neighbours. The station- 
ary state is trivial and non-interacting, while the dynamics is 
highly correlated and glassy when k» SI, 7. 



are understanding the interplay between static and dy- 
namic properties, the relevance of quantum versus clas- 
sical fluctuations, and the emergence of spatial correla- 
tions in the relaxation dynamics. In this work we take a 
first step towards addressing these issues from a dynam- 
ical facilitation perspective by introducing and studying 
a class of open quantum lattice systems that are the 
quantum analogs of classical facilitated spin models of 
glasses. These quantum spin systems, which are free of 
quenched disorder and have local interactions, display 
complex glassy dynamics despite the fact that their static 
properties are trivial. They also show a quantum version 
of dynamic heterogeneity, a feature we expect to be cen- 
tral to the dynamics of quantum glasses in general. 
Quantum facilitated models. Our aim is to construct 
a class of strongly interacting dissipative quantum many- 
body systems with the following two properties: firstly, 
interactions do not play an essential role in the statics 
and, as a consequence, the steady state of the system 
has the form of a direct product of single-site density 
matrices, like in a non-interacting problem. Secondly, the 
interactions make the dynamics and, in particular, the 
relaxation to the steady state, highly spatially correlated. 
Let us consider the lattice system of Figure [T] On each 
lattice site i — 1, . . . , N there is a quantum two-level sys- 
tem, and we identify the two levels with basis states |0)j 



and |l)j. These two levels are coupled coherently via 
an on-site Hamiltonian hi — £laf, where fi is the coher- 
ent coupling strength and erf = |0)j(l| + |l)i(0|. Each 
two-level system also interacts with a thermal bath that 
can cause incoherent excitation and de-excitation with 
rates 7 and k, respectively. In the absence of interaction 
between the sites, and under the standard Markovian ap- 
proximation for the thermal bath [13], the open quantum 
dynamics of each site is described by the master equation 

Pi = Ci{pi) = -i [h, pi] + na~ Pl o+ - - {afar, p^ 
+ja+ Pl ar - _ {a- a+ , Pi ) , (1) 

with a i being the ladder operators defined in terms of 
the usual Pauli matrices as af 1 = (erf ± iaf) /s/2. This 
equation has as stationary solution pf = \ u \u)i(u\ + 
A e |e)i(e|, where the unexcited (\u)i) and excited (|e)j) 
states can be written as a superposition of |0)j and |l)j 
and A Ui e are their corresponding stationary occupation 
probabilities [II]. When the two-level systems are non- 
interacting the stationary density matrix for the whole 
system is simply the direct product 
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We now convert the system into an interacting problem 
with the the same trivial stationary state. We introduce 
interactions between nearest neighbouring sites in the fol- 
lowing manner: we condition the rates for both coherent 
and incoherent changes at the i-th site to the state of 
neighboring sites, i.e., tt — > faQ, k —$■ fiK and 7 — > fij, 
where fi is a projection operator that acts only on the 
nearest neighbours of i, see Fig. [l] The many-body mas- 
ter equation now reads 



q = W( e ) 



with Hamiltonian and Lindblad 13 operators given by 
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The operators fi represent kinetic constraints. Under 
certain general conditions for /,, which we derive in the 
Supplemental Material [M] , the stationary state of the in- 
teracting problem is the trivial one (12J), i.e. g^ = g^ { . 
Here we focus on two specific choices for the constraints 
which define the two models we study in detail (in di- 
mension one, as higher dimensional generalisations are 
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FIG. 2: Kinetic constraints. For the classical models the up- 
per panel shows three sites in a configuration, assuming all 
other sites are in the state. For the quantum models the 
bottom panel shows three sites of a product basis state for 
the density matrix, where we use the notation u = \u)(u\ and 
e = |e)(e|. Red (dark on white) symbols represent sites that 
cannot evolve. Green (white on dark) symbols represent sites 
where the state of neighboring sites allow for their evolution. 
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qEast : 
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with P^ = |e)j(e| being projectors on the excited state 
on the i-th site. 

The lattice models defined by Eqs. ©-(pi) are quantum 
versions of the facilitated spin models for classical glasses 
7J. In particular, the choice of kinetic constraint (pi) 
defines a quantum Fredrickson-Andersen (qFA) model, 
while that of Q a quantum East (qEast) model [7] . The 
classical Fredrickson-Andersen (FA) and East models are 
recovered in the limit of vanishing coherent coupling, Q = 
0. Figure [2] sketches the effect of the kinetic constraints. 
For the classical FA and East models (upper panel), a 
site cannot change if both neighbours are in the state 0. 
In the FA model an excitation can facilitate changes to 
its left or its right, while in the East model facilitation is 
directional (e.g. only to the right as shown in the figure). 

In analogy to the classical case, in both the qFA and 
qEast models (bottom panel in Fig. [2]) there can be no 
evolution in a site surrounded by unexcited neighbours 
[see Eqs. (5J6)]. In the qFA an excitation, like the central 
one in the rightmost sketch, facilitates both its neigh- 
bours which evolve with the single site master operator 
of Eq. (fTl). In contrast, in the qEast the central excitation 
can only facilitate dynamics to its right. 

It is interesting to note that, in addition to the state 
of Eq. (pi), there is a second stationary state due to 
the kinetic constraints. This pure state, g( dark ) = 
\u ■ ■ ■ u)(u- ■ ■ u\, is a "dark state" which is annihilated 
by the Hamiltonian and all the Lindblad operators inde- 
pendently, since figl***) = £ (dark) /i = for all i. The 
dynamics defined by pH) is therefore reducible: It has 
one irreducible and "inactive" partition composed solely 
of g( dark ) and a second irreducible and "active" parti- 
tion composed of everything else. In the thermodynamic 
limit, N — > 00, the stationary state of the active par- 
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FIG. 3: Structure in trajectory space and dynamic heterogeneity. Quantum trajectories of the qFA, qEast and unconstrained 
models, from QJMC simulations for a ring of N = f spins. We plot (P^ )(t) S [0, 1] = [white, black] for all sites i — 1 to N. 
The rates (re, £7, 7) for the trajectories shown are: (1,0.25,0) for the qFA and (1,0.4,0) for the qEast and the unconstrained 
ones. Space-time correlations are evident in the trajectories of the constrained models, and absent in the unconstrained system. 



tition becomes g^ up to corrections which vanish ex- 
ponentially with N. In what follows we only consider 
dynamics in this partition. Nevertheless, the existence 
of the inactive state has important consequences for the 
dynamical behaviour of the system as we discuss below. 
Dynamical heterogeneity and slowdown of the re- 
laxation. We are interested in the dynamics of our quan- 
tum facilitated spin models in the regime where the co- 
herent coupling is weak and/or the temperature is low, 
i.e. k»J],7. In this regime X u ^> A e (see [2]), that is, 
the stationary state is such that very few sites are in the 
excited state |e). This however leads to a conflict with 
the dynamics, as the kinetic constraint on a site vanishes 
unless some of its neighbours are in the excited state, and 
in this regime most sites will be surrounded by unexcited 
sites. The consequence is a pronounced slowdown and 
the emergence of collective dynamics. 

Let us consider how excitations propagate in the quan- 
tum models we propose. Typically, at low CI/k, "//k, ex- 
citations will be isolated. We denote a state with such an 
excitation at site i by • • • u i _x e -i u i+\ ■ ■ ■ ( see Fig- [2])- As 
we have explained, the i-th site cannot evolve, but it can 
facilitate evolution of its neighbours. Let us first consider 
the case of the qFA model: here, the initial excitation 
can facilitate the excitation of either of the neighbouring 
sites. This new excitation subsequently facilitates the 
original site, which can now de-excite. The first process 
could be either virtual or real, so is limited by the rates 
ft or 7. The second process will be dominated by inco- 
herent de-excitation of rate k. In the classical limit, the 
outcome of such sequence is that an isolated excitation 
hops one site. In the qFA model excitations also propa- 
gate diffusively, with an effective diffusion constant that 
vanishes in the limit CI/k, j/k — > 0, in analogy to the 
classical FA model at low temperatures [7] . Dynamics of 
the qEast model is even more complex. Due to the direc- 
tionality of the kinetic constraint, Eq. (pi), the last step 
in the hopping sequence sketched above is not allowed 
and the relaxation is therefore hierarchical [2J [7] . 

Figure [3] shows examples of quantum trajectories from 
quantum jump Monte Carlo (QJMC) simulations [T5] , 
The trajectories for the qFA clearly show the diffusive 



nature of the dynamics: when initialized from a single 
site in the |e) state and the rest in the \u), this initial ex- 
citation propagates in a diffusive manner. Occasionally, 
excitations branch out, or coalesce. Due to the kinetic 
constraints excitations have to form connected chains in 
space and time [3J, as illustrated in the figure. The hier- 
archical relaxation dynamics of the qEast model is also 
observed in the trajectories of Fig. [3] For comparison 
we also show trajectories under sinrilar conditions for 
the unconstrained problem [see Eq. Q]. These trajec- 
tories are featureless, as there is no interaction between 
the sites, who then evolve independently. Note that all 
three systems shown in Fig.|3]possess the same stationary 
properties, despite the fact that their dynamics (and the 
approach to such stationary state) is obviously different. 
It is instructive to compare the quantum trajectories of 
Fig. [3] to those of the corresponding classical models, see 
e.g. Fig. 2 of Ref. [16]- 

The effective diffusion constant for propagation of ex- 
citations in the qFA model can be estimated analytically 
when ft, 7 -C K by carrying out an adiabatic elimination 
of fast degrees of freedom, or by numerical simulation 
of the evolution operator ^ in small lattices (see Sup- 
plemental Material [H]). These results suggest that the 
diffusion constant reads 
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6 k 



4ft 2 
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8ft 4 

K, 3 



(7) 



In the limit of zero temperature, 7 = 0, the diffusion 
constant is D — 8ft 4 /« 3 to leading order in ft. Relaxation 
of a site initially in an unexcited state will depend on 
how long it takes for the nearest excitation to diffuse 
to it. If the distance is I, then the time to relax will 
be t(1) ~ I 2 /D. Figure \n&) shows that this diffusive 
argument accounts very well for the relaxation time r rc i ax 
obtained from QJMC simulations in the qFA model. As 
expected, the relaxation slows down very markedly with 
the decrease of quantum fluctuations. The inset to Fig. 
|i[a) shows that the relaxation time in the qEast model 
grows with decreasing ft even faster than in the qFA as 
a consequence of the more restrictive kinetic constraints. 
From the trajectories of the qFA and qEast shown in 
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FIG. 4: Glassy slowing down and dynamical correlations, (a) 
Average relaxation time r re i ax as a function of 0,/k (at 7 = 0) 
in the qFA model. The dots are results from QJMC simula- 
tions. The line is the theoretical expectation for relaxation via 
activated diffusion of excitations: since the initial state has a 
single excitation, in a lattice of size N with periodic bound- 
ary conditions we expect r rc i ax w (2/iV) X//=i r (0> where 
t(1) ~ l 2 /D. Inset: r re i ax for the qEast model, (b) Fluc- 
tuations in waiting times between quantum jump events for 
the same qFA system, quantified via the Mandel parameter 
Qivtandei- A value of Qiviandci > is indicative of correlated 
fluctuations in the dynamics. 



Fig. [3] we see that, in the constrained cases considered, 
when the dynamics is slow it is also heterogeneous. Re- 
gions of space empty of excitations are slow to relax, as 
they need excitations external to the region to propagate 
into it for the dynamics to take place. The larger these 
empty regions the longer they take to evolve. These slow 
"space-time bubbles" |16j make the dynamics fluctuation 
dominated. Their origin can be traced back to the dark 
state g( dark ) : even if this is never accessed from the ac- 
tive partition, at low Q/k, j/k mesoscopic regions devoid 
of excitations look locally like the dark state; it is these 
spatial rare regions that give rise to dynamic heterogene- 
ity. Furthermore, in the limit of N — > oo, dynamic bot- 
tlenecks can become so pronounced that these systems 
display non-equilibrium (or trajectory) phase transitions 
between active and inactive dynamical phases, of the kind 
discussed in [17] . 

A manifestation of dynamic heterogeneity is in the 
fluctuations of the waiting times t w between quantum 
jump events in each site. When the dynamics is spa- 
tially correlated, quantum jumps are non-Poissonian and 
the corresponding distribution of waiting times is non- 
exponential. This can be quantified by the Mandel-Q 



parameter for the waiting times, Q 



Mandel 
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[18] . For correlated (bunched) quantum jumps we have 
that QMandci > 0, indicating that the fluctuations in the 
waiting times can be much larger than what is expected 
from a Poisson process. Figure W[b) shows that this is 
the case for the qFA model. In the glass literature this 
is often referred to as "decoupling" [2J [5] . 
Interplay of quantum and classical fluctuations. 
For the system to relax it is necessary to excite sites, 
which can be done incoherently, i.e., classically, with rate 
7, or coherently with rate f2. It would then seem that 
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FIG. 5: Enhanced glassiness due to weak quantum fluctua- 
tions. Effective diffusion constant D in the qFA model at 
7 7^ 0: It initially decreases as the strength of quantum fluc- 
tuations is increased from the classical limit (Do at fi = 0). 
Weak quantum fluctuations hinder, rather than promote, re- 
laxation. Inset: Similar behaviour of the branching rate B. 



adding quantum fluctuations to a classical facilitated sys- 
tem should always aid relaxation. But this is not the 
case. From Eq. (J7| we see that in the qFA model the 
first correction to the classical diffusion constant is nega- 
tive, and only for larger J7 quantum contributions become 
positive; see Figure [5] This means that weak quantum 
fluctuations enhance dynamical arrest by decreasing the 
propagation rate of excitations. The inset to Fig.[3]shows 
that this also holds for the rate B for branching of ex- 
citations [14]. This phenomenon of enhanced glassiness 
due to weak quantum fluctuations near the classical limit 
appears to be similar to that reported recently for glass 
forming quantum liquids in Ref. 12 . 

Outlook. We expect the dynamics of glassy quantum 
many-body systems with local interactions to display the 
essential feature of the idealized models we introduced in 
this paper, namely a spatially fluctuating and heteroge- 
neous dynamics which cannot be simply inferred from 
static behaviour. To understand this correlated dynam- 
ics one would need to uncover the effective kinetic con- 
straints. In this sense, we envisage quantum facilitated 
models as an effective description of the more complex 
glassy many-body quantum dynamics. Given the recent 
progress in quantum many-body physics with cold atoms 
it may be possible to directly implement and study con- 
strained models in experiments. For example, Rydberg 
atoms [T9] exhibit strong interactions that have been 
shown to lead to constraints in the coherent dynamics 
[20H22] which, in conjunction with atomic decay, can give 
rise to pattern formation in quantum jump trajectories 
[23] . Furthermore, recent work focusing on the engineer- 
ing of system-bath coupling [241 US] outlines a path for 
devising the constrained many-body jump operators re- 
quired by the models presented here. 
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no. EP/I017828/1 and Leverhulme Trust grant no. 
F/00114/BG. B.O. also acknowledges funding by Fun- 
dacion Ramon Areces. 
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SUPPLEMENTAL MATERIAL 

Non-interacting problem. Here we obtain the sta- 
tionary density matrix p\ s corresponding to Eq. (hj). It 



reads 

J st ) _ 



1 



8ft 2 + ( K + 7 ) z 



4ft 2 +7(^ + 7) -2«fi(«-7) 
2if2(K-7) 4^ 2 + k(k + 7) 



and can be diagonalized such that p). s = \ u \u)i(u\ + 
A e |e)i(e|. For k ~S> 7, ft, and to leading order in ft/re, 7/K 
one obtains that 



A>, 



7 16ft 4 , 7 16ft 4 

1-- ^-, A„«- + ;-, 
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2in(n — 1) 



1 - 



2fr 



2£1(k- K 7) 



Note that, in the limit of ft — > 0, the two states \u) and 
|e) can be identified with |0) and |1), respectively. 

Kinetic constraints. We are interested in introducing 
a constraint in the model such that the stationary state 
solution is the same as in the unconstrained case. This 
can be done by means of the transformation 

N 
i=\ 

where /, represents the kinetic constraint. The new mas- 
ter equation can be written as 

JV 



= Y.if^Mfl + UQ) 



with 

Me) = 



= -in 



(fie?[flQ]+[fi>Q]°fti 
%({*tf!fi<r7,Q}-fi{<rt<rT,Q}fl 
l({<rrflf«rt,Q}-fi{<rr<rt,<t}fl 



For this master equation to posses the same station- 
ary solution as the unconstrained problem, we need that 
Ai(Qni ) = 0' which is accomplished only if the following 
three properties of the constraints are fulfilled: (i) they 
are projectors, /, = fj = ff or, more generally, normal 
operators fifj — fjfi] (ii) the constraint at site i com- 
mutes with all operators at that site, [/;, of] — [/i, of ] — 
[/i,of] = 0; (hi) all constraints commute with the non- 
interacting stationary density matrix, [fi, g^ ] = Vi. 
The constraints of Eqs. (5J6) which define the qFA and 
qEast models satisfy these properties. 

Condition (iii) is not present in classical models; it is 
required as quantum density matrices in general are non- 
diagonal. It strongly restricts the class of operators that 
can be taken as kinetic constraints. In particular, on 



each site they have to be diagonal in the {\u), |e)} basis. 
Equations (3]|4) and conditions (i-iii) define the class of 
quantum facilitated spin models, of which the qFA and 
qEast are two representatives. 

Diffusion in the qFA model. Here we will obtain 
an analytic approximation for the diffusion constant in 
the qFA model in the regime where k 3> ft, 7 using the 
method of the adiabatic elimination of fast variables [55] . 
Let us consider the case of a system formed only by 
two sites, so that the basis is {\uu) , \eu) , \ue) , |ee)}. As 
discussed in the main text, the inactive state \uu) is dy- 
namically disconnected from any of the other states of 
the basis. Let us assume that initially no population is 
in this state, so that one can look at the dynamics of 
the reduced density matrix spanned only by the states 
{|ee) , \eu) , |ue)}. The master equation of the system is 
given by ([3]), and we divide it into two parts, 

g=[W K + W ](g), 

with the first term being proportional to the decay rate 
n which we assume is much larger than 7, ft and the 
second term being, at most, proportional to 7 and ft. 
In particular, one can see that the action of the first 
operator on the density matrix is given by 



.(<?) = 



We now proceed with the adiabatic elimination of fast 
variables, i.e., we look at the time evolution of the sys- 
temat typical time scales much larger than k _1 . First, 
we introduce the projection operator V = limt_>. 0o e Wre * 
that projects on the null space of W K , such that a = Vg 
is the stationary solution of a — W K (a) = 0. From the 
original master equation, and assuming that initially the 
state |ee) is not populated, the problem can be reduced 
to the effective time-evolution of the 2x2 density matrix 
a expressed in the basis {\eu) , \ue)}. The time-evolution 
of a is given by a = W Q a with 

W Q ='PWo-7'Wo[QWo + W K ] _1 QWo, (A.8) 

where Q = 1 — V. We can obtain from here, at least 
formally, the exact solution of the N = 2 problem which, 
in the the large k limit yields 

Q » [0i + 2 + 3 + 04] a 

with the operators Ok representing the contribution of 
order k in ft and 7: 0i = TW , 2 = A, 3 = B and 
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Up to fourth order in ft and second in 7, the effective 
equation of motion is given by a = W e ffa with the effec- 
tive master operator written in matrix form as 
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Rewriting this equation equivalently in terms of the 
site occupied by the excitation, i.e. |1) = \eu) and |2) = 
\ue) yields 

2 
a = ^D i aD\--i [ D\D h o^+R i aR\ 

- i{4^,a}+AaLl-i{xt£,, a }. 



Here, the operator Di = yPf[\i) {i\ [with 7i = 2" + 

^k 7§" (^^ 2 +7 2 )] corresponds to the dephasing of 

the two states of the spin on site i, and the operators 

Ri = i/DJ2\i) (i + 1| and L i = VW2 \i) (< ~ M to the 
hopping of the excitation to the right and left, respec- 
tively. The time-evolution of the population of each 
site, i.e., the diagonal elements of the density matrix 
a 1 = (i\a\i), is independent from the off-diagonal ones 
and read 



d t u l 



V 



[a l+1 



2a l 



In the continuum limit (a* — > a(x)), we obtain the fol- 
lowing Fokker-Planck equation 

da(x) d 2 a(x) 



Of 



V- 



dx 2 



and we can identify the parameter I? as a diffusion con- 
stant. 

Note, that we have implicitly assumed that the case 
of two atoms is sufficient to explain the diffusion for a 
large number of atoms. In practice, there are further 
corrections to be included due to the finite size of the 
system. When performing a similar analysis for larger 
systems (up to N = 4 sites), one finds that the actual 
form of the diffusion constant (up to fourth order in ft 
and second in 7) is given by the Eq. fcfj. 

Finally, in order to test the validity of the analytic re- 
sults, we have calculated numerically the master operator 
(A.8 1 for small systems (up to N = 5 sites) and obtained 



then the diffusion constant (D) and the branching ratio 



(B) as the element of W Q connecting euu ■ • ■ —} ueu . . . 
and euu •■■—>■ eue . . . , respectively. The results of these 
simulations for the diffusion constant and their compar- 
ison with the analytic solution are shown in Figure [6j 
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FIG. 6: The upper panel shows the diffusion constant from the 
numerically obtained master operator W Q l |A.8| l for lattices 
of size N = 2,3,4,5. The bottom panel shows collapse of 
the curves when the classical contribution is subtracted, D — 
Dq(N), where Do is the value of the diffusion constant at fl = 
0. The values used are: £> (2) = 7/2, D (3) = 7/2 + 7 2 /12k 
and A, (4) = D (5) = 7/2 - 7 2 /6k. 



